Forming giant planets via fragmentation of protoplanetary disks 
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ABSTRACT 



The evolution of gravitationally unstable protoplanetary gaseous disks has 
been studied using three dimensional smoothed particle hydrodynamics (SPH) 
simulations with unprecedented resolution. We have considered disks with 
initial masses and temperature profiles consistent with those inferred for the 
protosolar nebula and for other protoplanetary disks. We show that long-lasting, 
self-gravitating protoplanets arise after a few disk orbital times if cooling is efficient 
enough to maintain the temperature close to 50 K. The resulting bodies have 
masses and orbital eccentricities remarkably similar to those of observed extrasolar 
planets. 
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1. 

About 100 extrasolar planets have been detected by the wobble they induce on their star(i, 2). Their 
masses range from about one Jupiter mass (Mj) to more than 10 Mj and have orbits ranging from nearly 
circular to very eccentric. In the standard core-accretion model giant planets might require longer than 
10 6 years to form(3,^), which could exceed observed disk lifetimes (5 — 7). In particular, more than 80% 
of the stars in the Galaxy probably formed in dense clusters like those in the Orion nebula(#) where the 
ultraviolet radiation of bright stars can ablate the gaseous disk in far less than a million years(5, 6). Hence 
giant planet formation must occur quickly or such planets would be rare. Even in the case where a large 
solid core is assembled rapidly enough, torques acting between the disk and the protoplanets are believed 
to induce its complete inward migration in a few thousand years(P, 10) — planets could sink towards the 
star before being able to accrete the large gaseous masses observed (11, 12). Alternatively, giant planets 
could coagulate directly in the gas component as a result of gravitational instabilities in a cold disk with 
a mass comparable to that adopted in the core-accretion vno&e\(13 , 14)- Simulations done with codes that 
solve the hydrodynamical equations on a fixed grid show that slightly perturbed disks form strong spiral 
arms and overdensities at R > 10 AU(i5, 16) where the temperature can be lower than 60 K(i7, 18). The 
trigger of the instability might come from material of the protostellar cloud infalling onto the disk(13). If 
these condensations are long-lasting and can contract to planetary densities, gravitational instability would 
be the prevailing formation mechanism for giant planets because it takes less than a thousand years(i5, 15). 
Solid cores with masses as low as currently estimated for Jupiter (between and 10 Earth masses(iP)) 
could then form inside the gaseous protoplanets due to dust and planetesimals driven there by local pressure 
gradients in a few thousand years(20). 

However, due to the limitations of the techniques, simulations have not yet been able to show 
convincingly that the overdensities are not sheared apart by the tidal field of the star, nor that they can 
collapse into protoplanets(i6). One needs to achieve a high spatial resolution for such a purpose. Smoothed 
particle hydrodynamics (SPH) simulations ( 21 ) describe the gaseous medium as a collection of particles and 
can follow arbitrarily high densities. Disks have been simulated with this technique in the past but with 
less than 10 5 particles (22 — 24); with such a low mass resolution the evolution of the density distribution 
is considerably noisy and artificial fragmentation can take place(l?3). 

Here we report on the results of new 3D SPH simulations of marginally unstable disks of molecular 
hydrogen using as many as 1 million particles. The disks extend from 4 to 20 AU initially and arc in nearly 
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keplerian rotation around a solar mass star represented by a point mass. They have a minimum Toomre 
Q parameter of either 1.4 or 1.75, and masses of, respectively, 0.1 and 0.08M Q (see caption of Fig.l for 
details). In the initial stage the disks are evolved using a locally isothermal equation of state; the initial 
temperature decreases with radius following a power law profile predicted by detailed calculations of thermal 
balance between the central star, the disk and the protostellar cloud(i7) and is then held fixed locally. 
This approximation is based on the assumption that the cooling time is so short that the disk radiates away 
any thermal energy injection on a time-scale shorter than the orbital time. This seems to be supported by 
recent grid-based simulations that include radiative transfer in the diffusion approximation^ 6) but neglect 
the irreversible heating that can be generated by shock waves in a strongly unstable disk(^5, 26). While 
further investigation on the balance between heating and cooling in realistic disks will be needed in the 
future {26), here we concentrate on showing that actual protoplanets can form if the disk remains cold long 
enough(^7). 

After ~ 150 years, corresponding to about 5 orbital times at a radius of 10 AU, the disks develop 
trailing spiral arms and local overdensities at R > 10 AU. In the lighter disk the spiral arms grow in 
amplitude up to about 300 years, and then they settle down to a nearly stationary pattern (Fig.l). In the 
most massive disk (Q = 1.4) a two-armed mode grows in amplitude up to the point where, after about 
200 years, fragmentation occurs along the arms, and more than one distinct clump appears (Fig.l); then 
additional strong arms appear even at R ^ 7 — 8 AU, and more clumps are formed. Inside this radius the 
disk is too hot for condensations to form. Clumps quickly contract, reaching central densities of more than 
10 5 times the local density in a matter of a few orbital times (tens of years). These condensations are 
self-gravitating, their masses being larger than the local Jeans mass (28), and easily resist stellar tides. 

Such dense objects would be optically thick and would be unable to cool radiatively as efficiently as 
assumed by the locally isothermal approximation^ 6, 25). Therefore, we run again the same initial conditions 
changing the equation of state to adiabatic(£9) as soon as the spiral modes approach fragmentation, namely 
when they reach a density around ten times higher than the initial local density(-Z6, 26). The equation of 
state is changed throughout the disk. After 350 years non-axisymmetric features in the disk are weaker 
compared to the isothermal simulation, yet clump formation has proceeded. Nearly as many clumps as 
in the locally isothermal run are still present at R > 10 AU and have central densities still 10 5 times 
higher than the local density. The clumps are rapidly rotating spheroids; assuming conservation of angular 
momentum, bodies of about a Jupiter mass would have a rotation period of a few hours if they were allowed 
to contract further and reach the density of Jupiter. 
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The masses of the clumps shortly after all of them are in place (after ~ 350 years) range from 1 to 
5 Mj (all gas at densities at least ten times higher than the initial local density is identified as a clump), 
comparable to the masses of extrasolar planets (./). Clumps formed along the same spiral arm collide and 
merge into a more massive object sooner or later. In addition, all clumps accrete gas from the disk and 
begin to clear gaps. 

Time-steps can be as small as a few hours inside the clumps, and this slows down the simulation 
considerably. To keep following the clumps on a longer timescale we resort to simulations with a resolution 
5 times lower. By 350 years 7 clumps with masses comparable to the largest among the 13 clumps present 
in the higher resolution simulation have formed and are then followed for about 1000 years. We run both a 
locally isothermal and an adiabatic simulation to determine the role of the thermal structure of the disk in 
shaping the dynamics of clumps. 

Clumps are born on orbits with a wide range of eccentricities driven by the underlying strongly 
non-axisymmctric disk potential. These orbits then evolve in a variety of ways (Fig. 2). The smoother disk in 
the adiabatic run leads to less eccentric orbits; also, the higher pressure of the gas reduces considerably gas 
accretion by clumps. Slightly more circular orbits reduce the rate of close encounters and mergers between 
clumps relative to the isothermal run. After nearly a thousand years three clumps are left in the adiabatic 
run (Fig. 2) as opposed to two clumps in the isothermal run. These numbers are comparable to those of 
extrasolar giant planets in multiple systems (1). The masses of the surviving clumps, located between 3 
and 20 AU, are in the range 2-6 Jupiter masses in the adiabatic case and twice as big in the isothermal 
case. Most of the clumps have orbits with final eccentricities between 0.1 and 0.3 (Fig. 2) as many of the 
observed extrasolar planets (Fig. 2). The simulated protoplanets enter a rather quiescent evolutionary phase 
(no more mergers occur) several orbital times before the end of the simulations. However, inward migration 
might continue on time-scales longer than those explored heve(10 — 12), and orbits can also change due to 
interactions between the planets (30). 

This work shows that gravitational instability can actually form self-gravitating protoplanets and that 
long-lived systems with masses and orbits consistent with those of extrasolar planets arise. All this requires 
is to start with a marginally unstable disk (Q m in = 1.4 — 1.5) in which radiative cooling is efficient during 
the initial growth of the overdensities. We tested that hotter disks starting from a considerably higher 
Qmin (~ 2) can also become strongly unstable and form clumps if they are slowly cooled to temperatures 
comparable to those used in the disks that start with Q m i n ~ 1A{31). Therefore, clump formation does 
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not depend on how the disk reaches the state used in our initial conditions. Ice giant planets, like Uranus 
and Neptune, might also be formed by the same mechanism after a strong ultraviolet flux from nearby 
bright stars has photoionized the envelopes of protoplanets more massive than Jupiter, leaving a mostly 
metallic core(20). Because instabilities occur quickly, future observations of planets around very young 
stars will be a test for this model. In addition, direct imaging of giant planets at large distances (R > 
50 AU) from the stars(3,2) could also provide support to this model; as the outermost regions of the disk 
would be even cooler (and surface density and angular velocity fall equally with radius, i.e. as r~ 3 / 2 ) Q 
will be still decreasing at R > 20 AU, and thus fragmentation should occur out to these large distances. 
On the contrary, Jupiter-like planets would not form at such large distances in the core-accretion model 
because coagulation of planetesimals into a solid core would take too long with the small surface densities 
involved (33). Future observations of the gaseous medium in disks at different evolutionary stages, for 
example with the Space Infrared Telescope Facility (SIRTF), will help constrain the evolution of disk 
structure and will show whether it is consistent with the gravitational instability picture. 
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Fig. 1. — Snapshots of the simulations showing the protoplanetary disks seen face-on at different 
times. The colour-coded density on a logarithmic scale is shown out to 20 AU. Brighter colors trace 
higher densities; the density ranges between 10 -14 and 10 -6 g/cm 3 is shown using a logarithmic 
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scale. The evolution of two disks having a different initial minimum Toomre Q parameter, Q m i n , 
are shown. The Toomre parameter at a given disk location is defined as Q = Qv s /ttGT,, where E 
is the gas surface density, 0, is the angular velocity, G is the gravitational constant, and v s is the 
sound speed, v s = y/P/p, where P is the pressure and p is the density of the gas. Disks have a 
surface density profile E ~ r~ 3 / 2 (12) and temperature profiles as in Boss(9 — 11)- Q reaches its 
minimum at R > 10 AU, where the temperature is as low as 50 K, while Q > 4 close to the inner 
disk boundary, where the temperature is around 650 K. The two upper panels show the disk with 
initial Q min ~ 1.75 (M disk = 0.08M o ) at T = 160 yr (left) and T = 350 yr (right) , while the 
two lower panels show the disk with Qmin ~ 1.4 (Mdisk = O.1M ) at T = 160 yr and T = 350 yr. 
Gravity is softened on scales of 0.06 AU for disk particles. The central stellar potential is exactly 
keplerian at 2.5 AU and it is softened on smaller scales to speed up the computation. Both the 
central star and the inner disk boundary are free to move. The simulations were performed with 
GASOLINE, a parallel N-body/SPH code in which gravity is computed using a binary tree(3^). 
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Fig. 2. — Late stage of disk evolution. On top a face-on view of the disk after 800 years is shown 
for the 200.000 particles simulation in which the equation of state is switched to adiabatic after 
about 200 yr. The colour-coded logarithmic density (see Fig.l) is shown out to 25 AU. Three giant 
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protoplanets are left and their orbital evolution is shown in the plot at the bottom (thick lines) 
together with that of a clump from the simulation where the equation of state is kept isothermal 
(thin line). Each of the remaining clumps is the end result of a series of mergers; the orbital 
evolution of the most massive progenitor is shown. The complex combination of torquing by the 
non-axisymmetric disk and interactions with other clumps changes the orbital eccentricity and 
mean radius of the orbits. Overall the orbital evolution is considerably more complex than the 
nearly-steady inward migration expected in light, axisymmetric disks, eventually halting once the 
planet has cleared a gap (9 — 12). 



